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1  Introduction 


Background 

Integrated  Training  Area  Management  (ITAM)  Program 

The  Department  of  Defense  (DoD)  is  responsible  for  administering  more  than  25 
million  acres  of  federally  owned  land  in  the  United  States  (Public  Land  Law  Re¬ 
view  Commission  1970),  making  it  the  fifth  largest  Federal  land  management 
agency.  A  major  objective  of  the  ITAM  program  —  the  Army’s  program  for  man¬ 
aging  training  land  —  has  been  to  develop  a  method  for  estimating  training  land 
carrying  capacity.  The  Office  of  the  Deputy  Chief  of  Staff  for  Operations  and 
Plans  (ODCSOPS)  defines  training  land  carrying  capacity  as  the  amount  of 
training  that  a  given  parcel  of  land  can  accommodate  in  a  sustainable  manner, 
based  on  a  balance  of  use,  condition,  and  maintenance  practices.  ODCSOPS  has 
sponsored  the  Army  Training  and  Testing  Area  Carrying  Capacity  (ATTACC) 
program  to  estimate  this  training  land  carrying  capacity. 

ATTACC  Methodology 

The  ATTACC  methodology  is  used  to  estimate  training  and  testing  land  carrying 
capacity.  It  is  also  used  to  determine  land  rehabilitation  and  maintenance  costs 
associated  with  land-based  training  and  other  land  uses.  The  ATTACC  Hand¬ 
book  (U.S.  Army  Environmental  Center  [AEC]  1999),  Army  Regulation  (AR)  350- 
4,  and  Department  of  the  Army  Pamphlet  (DA  PAM)  350-4  document  the  stan¬ 
dard  operating  procedures  for  implementing  ATTACC. 

The  Evaluation  of  Land  Value  Study  (ELVS)  methodology,  a  precursor  to 
ATTACC,  was  an  initiative  sponsored  by  ODCSOPS  and  the  Assistant  Secretary 
of  the  Army  (Installations,  Logistics,  and  Environment)  [ASA  (IL&E)].  ELVS 
was  developed  to  estimate  training  land  carrying  capacity  and  the  cost  of  land 
rehabilitation  and  maintenance  associated  with  land-based  training  (Anderson 
et  al.  1996).  The  ELVS  methodology  quantified  training  land  condition  in  terms 
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of  Maneuver  Impact  Miles  based  on  mileage  projections  from  the  Battalion  Level 
Training  Model  and  training  event  and  vehicle  impact  severity  factors.  It  also 
used  a  modification  of  the  Revised  Universal  Soil  Loss  Equation  (RUSLE)*  to 
estimate  land  condition  in  terms  of  erosion  status  as  a  function  of  the  training 
load.  Land  rehabilitation  and  maintenance  costs  were  obtained  from  existing 
installation  records  and  regional  cost  estimates  of  particular  practices.  The 
ELVS  methodology  was  applied  to  eight  pacing  units  in  heavy  maneuver  train¬ 
ing  at  Fort  Hood,  TX,  and  the  Combat  Maneuver  Training  Center,  Hohenfels, 
Germany. 

The  ELVS  methodology  was  expanded,  updated,  and  redesignated  as  ATTACC. 
The  ATTACC  methodology  extends  the  ELVS  initiative  to  include  all  types  of 
Army  units  (including  unique  combat  units),  Army  service  schools,  Reserve 
Component  (RC)  units,  and  RC-unique  requirements.  The  ITAM  program  is  in¬ 
tegrating  the  ATTACC  methodology  into  the  Army’s  Weapon  System  Cost  Factor 
Development  Program  and  providing  tools  for  installation  personnel  to  use  to 
develop  local  requirements  and  impacts  analysis.  Portions  of  the  ATTACC 
methodology  were  incorporated  into  the  Range  Facilities  Management  and  Sup¬ 
port  System  to  provide  installation  schedulers  with  a  means  to  estimate  training 
load  during  scheduling  activities. 

ATT  ACC-Related  Army  User  Requirements 

Documentation  of  the  Army’s  environmental  technology  requirements  has  been 
an  iterative  process  that  began  with  a  series  of  meetings  in  1993  and  the  Office 
of  the  Directorate  of  Environmental  Programs’  (ODEP’s)  publication,  U.S.  Army 
Environmental  Requirements  and  Needs,  which  describes  the  critical  research, 
development,  test,  and  evaluation  needs  for  accomplishing  the  Army’s  mission 
with  the  least  impact  or  threat  to  the  environment.  These  Army-level  require¬ 
ments  were  reviewed  for  their  effect  on  readiness,  quality  of  life,  the  environ¬ 
ment,  and  the  timeliness  needed  for  the  Army  to  maintain  compliance  with  envi¬ 
ronmental  regulations.  All  major  commands,  major  subcommands,  ODCSOPS, 
and  the  Office  of  the  Deputy  Chief  of  Staff  for  Logistics  were  involved  in  estab¬ 
lishing  the  prioritized  and  validated  list  of  the  Army’s  environmental  technology 
requirements. 


Described  in  Chapter  2. 
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Land  Capacity  and  Characterization  is  the  third  priority  conservation  user  re¬ 
quirement  on  the  Army’s  list.  This  user  requirement  defines  the  Army’s  need  to 
estimate  training  land  carrying  capacity,  and  describes  the  ATTACC  methodol¬ 
ogy  as  a  means  to  provide  land  managers  with  scientifically  based  information  to 
support  sound  decisionmaking.  However,  this  user  requirement  also  defines  the 
current  version  of  ATTACC  as  limited  in  its  ability  to  provide  the  most  accurate 
information  for  decisionmaking.  This  limitation  is  related  to  the  accuracy  of  in¬ 
put  data  and  a  simplistic  characterization  of  the  three  components  of  the  model. 
The  user  requirement  identifies  research  and  development  required  to  improve 
the  accuracy  of  ATTACC. 

The  Land  Capacity  and  Characterization  user  requirement  identified  28  exit  cri¬ 
teria.  Each  criterion  defines  a  specific  product  required  to  address  a  specific  as¬ 
pect  of  the  overall  requirement.  Two  exit  criteria  address  the  temporal  and  spa¬ 
tial  uncertainty  of  model  inputs  and  predictions.  These  exit  criteria  are: 

1.  Develop  a  protocol,  tool(s),  and/or  factors  for  installation-level  use  that  reflect  a 
probable  range  of  results  in  the  ATTACC  methodology. 

2.  Develop  a  protocol,  tool(s),  and/or  factors  for  installation-level  use  that  improve 
spatial  results  from  ATTACC. 

ATTACC  Sensitivity  Analysis 

A  sensitivity  analysis  evaluates  the  magnitude  of  changes  in  a  model’s  output  as 
a  function  of  changes  in  the  input  parameter  values.  Moreover,  a  sensitivity 
analysis  of  a  model’s  responses  to  variations  in  input  values  can  be  used  to  de¬ 
termine  the  relative  importance  of  individual  input  values.  Results  of  a  sensitiv¬ 
ity  analysis  are  used  to  prioritize  data  acquisition  and  model  development  ef¬ 
forts. 

A  sensitivity  analysis  of  ATTACC  methodology  (Anderson  1999)  has  found  the 
methodology  to  be  sensitive  to  changes  in  the  climatic  factor  of  the  RUSLE.  This 
finding  implies  that  improvements  to  the  R  factor  portion  of  the  methodology  can 
result  in  an  overall  improvement  in  model  accuracy.  The  R  factor  will  be  de¬ 
scribed  in  the  next  chapter. 


Objectives 

The  purpose  of  this  research  was  to  determine  a  mapping  strategy  for  spatial 
and  temporal  prediction  and  uncertainty  analysis  of  the  rainfall  and  runoff  ero- 
sivity  R  factor.  The  strategy,  presented  in  this  report,  is  based  on  sequential 


10 


ERDC/CERL  TR-01  -39 


Gaussian  simulation.  Its  key  is  modeling  spatial  and  temporal  variability  of  the 
R  factor  using  semi-variograms  and  deriving  expected  estimates  and  their  vari¬ 
ances  at  any  location  and  time.  The  temporal  variability  of  this  factor  is  investi¬ 
gated  by  annual,  seasonal,  and  half-month  time  series  simulation.  The  specific 
objective  was  to  compare  uncertainty  using  the  simulation  and  traditional  iso- 
erodent  maps.  The  strategy  is  illustrated  in  a  case  study  in  which  no  rainfall 
data  were  available  for  the  area  of  interest.  The  study  area  was  expanded  to  in¬ 
clude  a  larger  area  where  rainfall  data  were  available  to  develop  prediction  mod¬ 
els.  From  the  larger  area,  data  for  the  area  of  interest  on  Fort  Hood  was  then 
extracted. 


Approach 

The  temporal  and  spatial  variability  of  the  RUSLE  R  factor  was  first  modeled 
using  semi-variograms  in  geostatistics.  Sequential  Gaussian  simulation  was 
then  used  for  uncertainty  analysis.  A  case  study  illustrates  the  strategy. 

Previous  applications  traditionally  used  R  factor  values  obtained  from  empirical 
isoerodent  maps  that  were  assumed  to  be  constant  over  time  and  space.  R  factor 
values  using  the  proposed  strategy  were  compared  to  the  traditional  approach  to 
demonstrate  differences  between  the  two  approaches. 


Scope 

The  project  presented  in  this  report  is  part  of  a  larger  initiative  titled  Error  and 
Uncertainty  Analysis  for  Ecological  Modeling  and  Simulation.  The  objectives  of 
the  larger  effort  are  to  develop  a  general  methodology  and  framework  for  spatial 
and  temporal  modeling,  simulation,  and  uncertainty  analysis  of  natural  re¬ 
source,  ecological,  and  environmental  systems.  Furthermore,  various  sources  of 
errors  such  as  measurement  errors,  sampling  errors,  model  errors,  expert  knowl¬ 
edge,  uncertainty,  etc.,  are  being  identified,  their  propagation  modeled  and  error 
budgets  developed  in  order  to  provide  management  guidelines.  The  case  study 
selected  for  this  overall  effort  is  the  ATTACC  methodology.  The  results  pre¬ 
sented  in  this  report  are  for  one  factor  of  the  RUSLE,  a  central  component  of 
ATTACC.  Similar  efforts  have  been  initiated  for  other  components  of  the 
RUSLE  and  ATTACC  models.  Results  of  these  efforts  are  being  reported  else¬ 
where.  The  information  provided  here  refers  to  the  ATTACC  methodology  as  de¬ 
scribed  in  the  ATTACC  Handbook  (AEC  1999). 
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Mode  of  Technology  Transfer 

This  strategy  for  RUSLE  R  factor  prediction  and  uncertainty  analysis  will  be 
provided  directly  to  Army  personnel  responsible  for  ATTACC  implementation. 
The  information  is  also  provided  to  organizations  responsible  for  developing  and 
refining  the  ATTACC  methodology. 


Units  of  Weight  and  Measure 

Some  U.S.  standard  units  of  measure  are  used  in  this  report.  A  table  of  conver¬ 
sion  factors  for  Standard  International  (SI)  units  is  provided  below. 


12 


ERDC/CERL  TR-01  -39 


2  Approaches  for  Derivation  of  R  Factor 
Values 


The  RUSLE  is  widely  used  to  predict  average  annual  soil  loss  for  specific  loca¬ 
tions.  RUSLE  is  an  empirical  model  in  which  soil  loss  depends  on  the  factors  of 
rainfall  erosivity  (R),  soil  erodibility  (K),  slope  steepness  (S),  slope  length  (L), 
cover  management  (C),  and  support  practice  (P)  (Renard  et  al.  1997). 

Soil  erosion  is  greatly  influenced  by  the  intensity  and  duration  of  precipitation 
events  and  by  the  amount  and  rate  of  the  resulting  runoff.  The  R  factor  is  a 
quantitative  expression  of  the  erosivity  of  local  average  annual  precipitation  and 
runoff.  R  factor  values  are  derived  from  research  data  (Wischmeier  1959; 
Wischmeier  and  Smith  1958).  Erosivity  index  (El)  values  are  estimated  for  each 
precipitation  event.  The  annual  rainfall  and  runoff  erosivity  factor  is  the  sum  of 
the  El  values  for  all  rain  events  in  a  year.  Larger  rainfall  and  runoff  erosivity 
factor  values  indicate  climatic  conditions  with  higher  potential  for  soil  loss. 

Isoerodent  maps  developed  using  historical  data  (frequently  a  30-year  period) 
are  widely  used  to  estimate  a  rainfall  El  for  a  specific  area  (Renard  et  al.  1997). 
Linear  interpolation  is  used  to  estimate  El  between  the  contour  lines  of  the  iso¬ 
erodent  maps.  This  approach  assumes  that  the  rainfall  erosivity  factor  is  linear 
over  space  and  constant  over  time.  As  suggested  by  McGregor,  Mutchler,  and 
Bowie  (1980),  this  assumption  may  not  be  true  in  all  cases.  For  example,  global 
climate  change  might  result  in  a  change  in  the  frequency  of  annual  precipitation 
and  thus  in  a  change  of  rainfall  El. 

Although  spatially  explicit  El  can  be  derived  by  linear  interpolation,  a  constant 
value  for  an  area  is  usually  applied.  This  application  may  lead  to  spatial  predic¬ 
tion  smoothing  and  unaccounted  uncertainty  related  to  this  smoothing.  Fur¬ 
thermore,  uncertainty  of  the  rainfall  El  from  the  isoerodent  maps  is  usually  un¬ 
known.  Additionally,  calculation  of  R  factor  values  for  new  locations  and  areas  is 
laborious  and  requires  long-term  rainfall  intensity  data. 

An  alternate  method  is  to  use  a  simple  linear  regression  between  the  rainfall  El 
and  average  annual  precipitation  data  (Bollinne  et  al.  1980).  Mikhailova  et  al. 
(1997)  further  developed  a  regression  of  the  rainfall  El  as  a  function  of  average 
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annual  precipitation  and  elevation.  Elevation  was  found  to  be  highly  significant 
in  predicting  rainfall  EL  Goovaerts  (1999)  presented  a  geostatistical  method  for 
this  purpose,  where  a  digital  elevation  model  is  incorporated  into  mapping  an¬ 
nual  and  monthly  erosivity  values.  His  study  suggested  that  geostatistical 
methods  were  promising  for  monitoring  and  mapping  rainfall  and  runoff  erosiv¬ 
ity.  Geostatistical  methods  are  divided  into  kriging  and  simulation  algorithms, 
which  are  extensively  applied  in  natural  resource  and  environmental  sciences. 

Kriging  is  a  generalized  least-squares  algorithm  minimizing  the  magnitude  of 
the  error  variance  with  unbiased  estimates,  and  can  be  used  for  making  spatial 
prediction  at  unknown  locations  based  on  measurements  from  known  locations. 
Modeling  prediction  uncertainty  is  based  on  semi-variograms  measuring  the  spa¬ 
tial  variability  of  the  data.  Kriging  estimates  are  smoothed  and  are  best  in  local 
prediction;  however,  kriging  variances  depend  only  on  the  data  configuration  and 
not  on  the  actual  observed  data,  and  thus  do  not  adequately  reflect  uncertainty. 
Examples  of  kriging  include  the  interpolation  of  heavy  metal  concentrations  in  a 
contaminated  site  using  three  kriging  estimators  by  Juang  and  Lee  (1998)  and 
the  work  of  Rogowski  and  Wolf  (1994)  on  variability  of  soil  map  unit  delineation 
using  a  kriging  method. 

The  other  geostatistical  approach  uses  spatial  simulation  techniques  whereby 
the  conditional  distributions  are  developed  according  to  collected  data  sets. 
From  these  distributions  the  values  of  the  stochastic  variable  at  unknown  loca¬ 
tions  are  drawn  at  random.  Once  values  at  all  the  unknown  locations  are  simu¬ 
lated,  a  realization  of  the  stochastic  variable  is  developed.  After  many  realiza¬ 
tions,  the  set  of  alternative  realizations  provides  a  visual  and  quantitative 
measure  (actually  a  model)  of  spatial  uncertainty  (Wang  et  al.  2000;  Goovaerts 
1997).  The  expected  estimates  and  various  uncertainty  measures  such  as  condi¬ 
tional  variances  and  probability  maps  can  be  derived  from  the  realizations.  The 
error  variances  depend  not  only  on  data  configuration,  but  also  on  the  data 
themselves,  thus  implying  global  uncertainty.  The  most  commonly  used  tech¬ 
nique  is  sequential  Gaussian  simulation. 


14 


ERDC/CERL  TR-01  -39 


3  Study  Area,  Data  Set,  and  Rainfall 
R  Factor 

Fort  Hood,  the  study  area  for  this  research,  occupies  87,890  hectares  (ha)  (Na- 
kata  1987)  in  the  central  Texas  counties  of  Bell  and  Coryell.  Fort  Hood’s  climate 
is  characterized  by  long,  hot  summers  and  short,  mild  winters.  Average  tem¬ 
peratures  range  from  a  low  of  about  8  °C  in  January  to  a  high  of  29  °C  in  July. 
Average  annual  precipitation  is  81  cm.  Table  1  shows  average  monthly  and 
maximum  rainfall  values  for  Fort  Hood. 

Elevation  at  Fort  Hood  ranges  from  180  to  375  m  above  sea  level  with  90  percent 
of  the  area  below  260  m.  Most  slopes  are  in  the  2  to  5  percent  range,  with  slopes 
in  excess  of  45  percent  occurring  as  bluffs  along  the  flood  plain  and  as  the  sides 
of  slopes  on  the  hills.  Fort  Hood  lies  in  the  Cross  Timbers  and  Prairies  vegeta¬ 
tion  area  (Gould  1975),  which  is  normally  composed  of  oak  woodlands  with  grass 
undergrowth.  Soil  cover  is  generally  shallow  to  moderately  deep  and  clayey,  un¬ 
derlain  by  limestone  bedrock. 

Based  on  traditional  isoerodent  maps,  the  annual  R  factor  value  for  Fort  Hood  is 
a  constant  270  (Renard  et  al.  1997).  No  rainfall  observation  stations  were  within 
the  study  area,  however.  It  was  necessary,  therefore,  to  use  data  from  rainfall 
observation  stations  in  the  surrounding  area.  The  247  rainfall  stations  ulti¬ 
mately  used  were  located  in  Texas  and  the  surrounding  states  of  Arkansas,  Colo¬ 
rado,  Kansas,  Louisiana,  New  Mexico,  and  Oklahoma.  Figure  1  shows  the  geo¬ 
graphical  location  of  these  stations.  Out  of  the  247  stations,  29  were  randomly 
sampled  and  used  as  a  validation  data  set.  The  remaining  218  stations  were 
used  to  develop  the  spatial  models. 


Table  1.  Average  annual  climatic  data  for  Fort  Hood,  TX. 


Month 

January 

5.1 

5.8 

February 

5.6 

6.3 

March 

5.6 

5.1 

April 

9.4 

4.3 

May 

11.4 

22.6 

June 

8.4 

8.6 

July 

4.6 

4.1 

August 

7.1 

6.1 

September 

8.9 

12.2 

October 

9.7 

7.9 

November 

5.1 

5.6 

December 

4.3 

3.8 
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The  R  factor  values  for  rainfall  erosivity  were  calculated  for  each  rainfall  station 
by  a  method  developed  by  a  research  team  headed  by  Hollinger  *  This  team  is  in 
the  process  of  calculating  R  factor  values  on  a  half-month  basis  across  the  entire 
United  States.  An  annual  R  factor  value  is  the  sum  of  rainfall  El  for  all  rain 
events  within  a  year.  If  a  rainstorm  produces  less  than  0.5  in.  of  rain  during  a 
period  of  6  hours  (Wischmeier  and  Smith  1978),  the  El  of  the  rainstorm  is  ob¬ 
tained  by  multiplying  kinetic  energy  of  the  rain  with  30-minute  maximum  rain 
intensity.  Equation  1  was  used  to  calculate  the  energy  contained  in  the  volume 
of  rain: 

Equation  1.  Erosivity  equation. 

E  =  0.29x[l  -0.72xexp(-0.082xl)] 

where  E  =  the  kinetic  energy  in  the  shower  (MJ  ha1  mm-1) 

I  =  the  shower  intensity  in  mm  hr1. 

The  calculation  of  an  R  factor  value  requires  rainfall  data  over  a  long  term,  usu¬ 
ally  more  than  20  years.  In  this  study,  a  data  set  of  21  years  was  used.  The  El 
for  each  year  was  derived  and  the  mean  value  for  all  the  years  was  considered  to 
be  the  rainfall  erosivity  R  factor.  Seasonal  and  half-month  rainfall  erosivity  val¬ 
ues  were  computed  using  the  mean  values  of  El  over  the  21  years  for  each  half¬ 
month  and  season. 


Steven  Hollinger  is  a  Senior  Program  Scientist  with  the  Illinois  State  Water  Survey,  Champaign,  IL. 
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4  Methods  for  Spatial  and  Temporal 
Modeling 

Experimental  semi-variograms  measuring  the  spatial  and  temporal  variability  of 
rainfall  erosivity  were  first  calculated  and  modeled  for  the  large  area  (Fort  Hood, 
Texas,  and  the  surrounding  states).  Spatial  and  temporal  prediction  and  uncer¬ 
tainty  analysis  for  annual,  seasonal,  and  half-month  R  factors  were  then  done 
using  sequential  Gaussian  simulation.  The  prediction  and  variance  maps  of  the 
R  factor  values  were  derived  over  space  and  time.  The  prediction  and  uncer¬ 
tainty  maps  for  the  Fort  Hood  study  area  were  extracted  from  those  for  the  large 
area.  Simulations  were  tested  using  the  differences  between  estimates  and  ob¬ 
servations  from  the  validation  data.  The  results  were  compared  with  those  ob¬ 
tained  using  traditional  isoerodent  maps  in  terms  of  uncertainty  analysis. 

Because  the  rainfall  stations  were  not  systematically  located  in  the  larger  study 
area,  data  for  the  rainfall  erosivity  R  factor  were  de-clustered  to  obtain  a  reason¬ 
able  weight  for  each  station.  Further,  the  Gaussian  simulation  requires  normal 
distributions  of  data,  so  normal  score  transformation  of  the  original  data  was 
done  to  make  the  transformed  data  normally  distributed.  The  spatial  variability 
of  the  transformed  data  (instead  of  the  original  data)  was  then  modeled  using 
semi-variograms  and  used  in  the  simulations.  The  methods  are  briefly  described 
below. 


Semi-variogram 

A  semi-variogram  measures  the  spatial  variability  of  a  random  variable.  By 
sampling  a  random  variable  (z)  in  a  study  area,  observations  (n)  are  acquired 
(z(ua)  (a  =  1,  2,...,n)).  The  vector  of  spatial  coordinates  of  the  ath  individual  is  ua. 
The  experimental  semi-variogram  g(h)  is  generally  calculated  from  samples  us¬ 
ing  the  following  expression  (Krige  1966): 

Equation  2.  Semi-variogram  equation. 

1  N(h) 

7(h)  =  ——  £  (z(ua)~  z(ua  +  A))2 
N(h)  SS 
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where  N(h)  =  the  number  of  pairs  used 
h  =  vector  separating  two  values 

z(ua)  and  z(ua+h)  =  the  two  values  of  the  interest  variable  separated  by  a 
distance  of  h. 

The  magnitude  of  the  semi-variogram  quantifies  the  spatial  variability.  The  lar¬ 
ger  the  magnitude  of  the  semi-variogram,  the  higher  the  spatial  variability  and 
the  lower  the  spatial  correlation  in  the  data  set. 

Ideally,  the  value  of  the  semi-variogram  should  be  zero  when  the  separation  vec¬ 
tor  h  is  zero.  In  practice,  this  is  usually  not  true  because  of  measurement  errors. 
In  this  case,  a  so-called  “nugget  effect”  (nugget  variance)  exists.  As  the  semi- 
variogram  increases,  separation  distance  increases  and  reaches  its  maximum  at 
a  distance  called  range  parameter.  The  maximum  semi-variogram  value  is  sill 
parameter  (structured  variance).  The  experimental  semi-variograms  are  often 
fit  using  several  models  that  include  spherical,  Gaussian,  and  exponential  mod¬ 
els,  and  the  model  that  best  fits  the  data  is  selected.  Additionally,  different  di¬ 
rections  are  taken  into  account  to  determine  whether  the  spatial  variability  is 
isotropic  or  anisotropic. 


Simple  Kriging 

Simple  kriging  allows  the  estimation  of  unknown  locations  based  on  information 
gathered  from  the  known  locations  within  a  neighborhood  given  a  radius  for 
searching  for  the  data.  Kriging  estimates  incorporate  the  spatial  variability 
structure  of  the  data  into  predictions.  The  simple  kriging  estimator  (Krige  1966) 
is  given  as: 

Equation  3.  Simple  kriging  estimator  equation. 


n(u) 

a- 1 

where  Z]k  (u)  =  a  simple  kriging  estimate  at  an  unknown  location  (u) 

n(u)  =  the  number  of  the  data  used  at  the  known  locations  given  a 
neighborhood 

A*  (u)  =  the  weight  assigned  to  the  datum  z(ua) 
m  =  the  expectation  value  of  Z(u). 


Mu) 


a= 1 
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In  simple  kriging,  the  mean  should  be  a  known  constant  over  all  the  area.  The 
error  variance  is  minimized  under  the  constraints  that  it  is  unbiased  with  an  er¬ 
ror  mean  of  zero.  The  error  variance,  <J2sk ,  is  usually  written  in  terms  of  covari¬ 
ance  functions,  C(.)  (Goovaerts  1997),  and  it  is: 

Equation  4.  Error  variance  equation. 


n(u) 

o1,i=C(0)-'^Kl{<‘)C(ua-u) 

a= 1 


Sequential  Gaussian  Simulation 

Sequential  Gaussian  simulation  assumes  that  the  underlying  distribution  is 
Gaussian.  As  such,  the  appropriateness  of  the  data  distribution  must  thus  be 
tested  before  simulation.  Distribution  that  is  not  normal  requires  a  normal  score 
transformation  of  the  original  data  into  a  new  data  set  with  a  standard  normal 
cumulative  density  function,  and  then  the  simulated  values  need  to  be  trans¬ 
formed  back  to  the  original  scale. 

A  data  set  {z(ua),a  =  l,2,3,...,n}  is  sampled  from  a  study  area  divided  into  N 
nodes  where  |z(Uj),j  =  l,2,3,...,N]is  a  set  of  random  variables  defined  at  N  loca¬ 
tions  {uj }.  Conditionally,  L  joint  realizations  |z(l)(Uj),j  =  l,...,N]  exist  for  1  =  1,  2, 

...,  L.  The  N  random  variables  can  be  generated  in  the  sequential  Gaussian 
simulation.  The  N-point  conditional  cumulative  density  function  is  expressed  as 
the  product  of  N  one-point  conditional  cumulative  density  functions  given  the  set 
of  n  original  data  values  and  N-l  realizations  (Goovaerts  1997).  If 
F(uN;zN  |  (n+N-1))  is  the  conditional  cumulative  density  function  of  Z(uN)  given 

the  set  of  n  original  data  values  and  the  N-l  previous  realizations 
Z(Uj)  =  z(l)(Uj),j  =  l,...,N-l,  the  N-point  conditional  cumulative  density  function 

follows: 

Equation  5.  Cumulative  density  function. 

F(uj,...,uN;z,,...,zN  |  (n))  =  F(un;zn  |(n+N-l))x 
F(un-i  5  zn-i  |(n+N-2))x  ...x 
F(uj ;  z2  |(n+ 1 ))  x  F(u;  z,  |(n)) 
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The  steps  for  the  simulation  are: 

1 .  Test  the  appropriateness  of  multi-Gaussian  distribution  by  declustering  and  per¬ 
forming  normal  score  transformation  of  the  original  data  and  developing 
standardization  semi-variograms. 

2.  Define  a  random  path  to  visit  each  node  of  the  grid  only  once  in  the  study  area. 

3.  At  the  first  node  to  be  visited,  calculate  the  mean  and  variance  of  the  Gaussian 
conditional  cumulative  density  function  given  the  n  original  data  using  simple 
kriging  and  the  modeled  normal  score  semi-variogram.  From  the  conditional  dis¬ 
tribution  above,  draw  a  value,  transform  it  back  to  the  original  data,  and  add  it  to 
the  data  set. 

4.  At  the  ith  node  to  be  visited,  determine  the  parameters  of  the  Gaussian  distribu¬ 
tion  given  the  n  original  data  and  all  (i-1)  simulated  values  at  the  locations  previ¬ 
ously  visited,  and  for  the  ith  node.  From  that  conditional  distribution  draw  a 
value  that  is  transformed  back  and  becomes  a  conditional  datum  for  all  subse¬ 
quent  drawings. 

5.  Repeat  step  4  until  all  N  nodes  are  visited  and  provided  with  simulated  values. 

Running  L  times,  each  time  with  a  possible  different  path  to  the  N  nodes,  will 
lead  to  L  realizations  from  which  an  expected  value  and  prediction  variance  for 
each  node  can  be  derived.  The  realizations  thus  provide  a  visual  measure  and  a 
model  of  spatial  uncertainty. 

As  previously  mentioned,  simple  kriging  estimates  are  smoothed  and  are  best  in 
local  prediction;  however,  kriging  variances  depend  only  on  the  data  configura¬ 
tion  and  not  on  the  actual  observed  data,  and  thus  do  not  adequately  reflect  un¬ 
certainty.  The  advantage  of  the  sequential  Gaussian  method  over  kriging  meth¬ 
ods  is  that  it  provides  spatial  uncertainty  analysis  by  calculating  conditional 
variances.  The  conditional  variance  depends  on  not  only  data  configuration  but 
also  on  data  values.  This  conditional  variance  in  theory  will  provide  a  more  real¬ 
istic  assessment  of  uncertainty  across  space. 
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5  Results 


Spatial  and  temporal  variation  of  annual  and  seasonal  rainfall  erosivity  R  factor 
observations  for  the  large  area  are  shown  in  Figure  2.  The  annual  and  seasonal 
R  factor  values  were  lower  at  the  west,  increased  from  the  west  to  the  central 
east,  and  slightly  decreased  further  to  the  east.  The  annual  R  factor  values  in¬ 
creased  slightly  from  the  north  to  the  central  and  then  decreased  slightly  to  the 
south.  Winter,  spring,  and  autumn  R  factor  values  increase  slightly  from  the 
north  to  the  south,  and  the  summer  values  decreased  slightly  along  the  direction 
of  the  north  to  south.  Overall,  the  summer  R  factor  values  were  the  largest,  fol¬ 
lowed  by  autumn,  spring,  and  winter  R  factor  values,  which  were  the  smallest. 
The  Fort  Hood  study  area  had  no  rainfall  stations;  however,  there  were  four  sta¬ 
tions  in  the  proximity  of  Fort  Hood.  Their  annual  rainfall  erosivity  R  factor  val¬ 
ues  varied  from  317  to  392. 

The  observations  of  half-month  rainfall  R  factor  values  are  not  presented  here 
because  of  space  limitations.  The  spatial  variability  of  the  observations  for  each 
half-month  map  was  similar  to  that  of  the  annual  and  seasonal  ones  in  Figure  2. 
That  is,  the  largest  values  were  distributed  at  the  central  east  and  the  smallest 
at  the  west,  but  their  spatial  variability  was  more  smooth  compared  with  that  of 
annual  and  seasonal  R  factor  values.  The  half-month  R  factor  values  and  their 
spatial  variability  rose  slowly  from  the  first  half  of  January  to  the  second  half  of 
March,  then  rose  rapidly  until  the  first  half  of  June.  After  that  period,  they  be¬ 
gan  a  trend  of  very  slight  decrease  and  fluctuation  until  the  first  half  of  October, 
then  decreased  until  December.  The  smallest  values  and  variability  occurred  in 
the  first  half  of  January  and  the  largest  occurred  in  the  first  half  of  June. 

Figure  3  shows  an  example  comparison  of  the  data  distribution  with  and  without 
normal  score  transformation.  The  distribution  of  original  data  was  obviously  not 
normal.  The  transformation  led  to  a  new  data  set  that  was  normally  distributed. 
The  spatial  variability  of  annual,  seasonal,  and  half-month  R  factor  values  were 
then  analyzed  by  experimental  semi-variograms  in  the  directions  of  azimuth  0, 
45,  90,  and  135  degrees.  The  semi-variograms  were  similar  at  all  the  directions 
and  considered  isotropic.  Figure  4  shows  an  example  of  isotropic  semi- 
variograms. 
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Annual  rainfall-runoff  erosivity 

Annual  R 

170  -  215 
215  -  260 
260  -  305 
'  305  -  350 

'  350  -  395 

»  395  -  440 
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Seasonal  rainfall-runoff  erosivity 

Winter  Spring 


l  • 


Summer 


Autumn 


Seasonal  R 

-  2  -  35 
35  -  68 
68  -  101 
101  -  134 
134  -  167 

-  167  -  200 
200  -  234 


Figure  2.  Spatial  location  of  rainfall  stations  in  the  large  area  used  for  calibration  and  their 
annual  (top  right)  and  seasonal  rainfall  erosivity  R  factor  (bottom  left). 


Figure  3.  Comparison  of  data  distribution  without  (left)  and  with  (right)  normal  score 
transformation  for  the  rainfall  erosivity  R  factor  for  the  second  half  of  September. 
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Figure  4.  An  illustration  of  the  examining  isotropic  feature  of  experimental  semi-variograms  for 
four  directions  (00,  450,  900,  and  1350)  for  the  rainfall  R  factor  during  the  second  half  of 


September. 


The  experimental  and  modeled  semi-variogram  of  annual  R  factor  values  using 
the  original  data  is  shown  in  Figure  5.  The  Gaussian  model  was  the  best  semi- 
variogram  model  in  terms  of  fit  for  the  original  data.  Figure  6  shows  experimen¬ 
tal  and  modeled  semi-variograms  of  seasonal  rainfall  R  factor  values  using  the 
original  data.  The  Gaussian  model  was  also  found  to  best  fit  the  semi- 
variograms.  With  the  normal  score  transformation  and  standardization,  the 
semi-variograms  for  annual  and  seasonal  rainfall  R  factor  values  were  further 
modeled.  Their  nugget,  sill,  and  range  parameters  are  in  Table  2,  together  with 
those  parameters  using  the  original  data.  After  standardization,  the  sum  of  sill 
and  nugget  was  forced  to  be  one  unit.  The  experimental  and  modeled  semi- 
variograms  for  half-month  rainfall  R  factor  values  were  developed.  The  figures 
and  model  parameters  are  omitted  here  because  of  space  limitation.  A  total  of  24 
semi-variogram  models  were  obtained,  including  17  Gaussian,  6  spherical,  and  1 
exponential  models. 


Figure  5.  Experimental  (dots)  and  modeled  (line)  omni-directional 
semi-variogram  for  annual  rainfall  erosivity  R  factor  values. 


Figure  6.  Experimental  (dots)  and  modeled  (line)  omni-directional  semi-variograms  for  seasonal 


rainfall  erosivity  R  factor  values. 
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Table  2.  Parameters  of  semi-variograms  using  original  and  standardized  data. 


Rainfall 

R  factor 

Original  data 
(nonstandardized) 

Normal  score  data 
(standardized) 

Sill 

B9 

Model 

KffiflU 

Sill 

Model 

Annual 

1085.6 

9233.6 

508.7 

G 

0.353 

0.647 

827.9 

G 

Winter 

106.6 

543.6 

558.5 

G 

0.150 

0.850 

900.0 

S 

310.8 

1011.9 

435.8 

G 

0.300 

0.700 

1000.0 

S 

Summer 

496.3 

965.1 

920.3 

G 

0.400 

0.600 

800.0 

S 

Autumn 

375.1 

775.9 

369.8 

G 

0.350 

0.650 

1200.0 

S 

Note:  G  =  Gaussian  model;  S  =  Spherical  model 


The  prediction  and  variance  maps  of  annual  rainfall  R  factor  values  for  the  large 
area  are  shown  in  Figure  7.  The  predicted  values  of  the  annual  R  factor  had 
similar  spatial  distribution  to  that  of  the  observations.  That  is,  they  increased 
from  the  west  to  the  east  and  from  the  north  to  south,  but  the  highest  values  ex¬ 
isted  in  the  central  east.  As  would  be  expected,  the  prediction  variances  were 
higher  outside  the  area  of  rainfall  stations.  In  the  large  area,  the  denser  the 
rainfall  stations,  the  lower  the  prediction  variance.  The  lowest  variances  existed 
in  the  central  east  area  because  of  the  high  density  of  rainfall  stations.  Further, 
variances  were  higher  in  the  central  west  due  to  the  lower  density  of  rainfall  sta¬ 
tions. 
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No  Data 


Figure  7.  Prediction  and  variance  images  of  annual  rainfall  erosivity  R  factor  values  for  the  large 
area  (left)  and  Fort  Hood  (right),  and  the  probability  map  for  annual  rainfall  erosivity  being  larger 
than  270  (bottom  right). 
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The  prediction  and  variance  maps  for  the  Fort  Hood  area  were  then  extracted 
from  the  corresponding  maps  for  the  large  area  and  are  shown  at  the  right  in 
Figure  6.  The  figure  also  shows  the  probability  map  for  the  annual  R  factor 
value  larger  than  270  at  Fort  Hood.  The  predicted  values  varied  from  350  to 
376.  The  estimates  were  not  a  constant  over  space,  however,  and  much  higher 
than  270  derived  from  the  commonly  used  isoerodent  map  for  Fort  Hood.  The 
probabilities  of  all  the  predicted  values  larger  than  270  were  higher  than  0.7. 
The  variances  of  the  predicted  values  ranged  from  3600  to  5437. 

Figure  8  shows  the  results  of  spatial  and  temporal  prediction  and  uncertainty 
analysis  for  seasonal  rainfall  erosivity  R  factor  values  for  both  the  large  area  and 
the  Fort  Hood  area.  The  predicted  values  in  the  large  area  (Figure  8,  upper  left) 
varied  over  space  and  time.  The  values  were  consistent  with  observations  of  the 
seasonal  R  factor  values  in  Figure  2.  The  lowest  prediction  was  in  the  west  and 
the  largest  in  the  central  east.  Winter,  spring,  and  autumn  prediction  values 
slightly  increased  from  the  north  to  the  central  and  slightly  decreased  from  the 
central  to  the  south.  For  the  summer,  high  prediction  values  were  in  some  local¬ 
ized  areas  in  the  north  and  south.  Over  the  four  seasons,  the  summer  had  the 
largest  predicted  values  and  the  winter  had  the  lowest.  The  predicted  values  in 
spring  and  autumn  were  similar. 
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Figure  8.  Prediction  (top)  and  variance  (bottom)  images  of  seasonal  rainfall  erosivity  R  factor 
values  for  the  large  area  (left)  and  Fort  Hood  area  (right). 
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The  variances  of  the  predicted  values  in  spring  and  autumn,  however,  were  lar¬ 
ger  than  those  in  winter  and  summer  (Figure  8,  lower  left).  The  winter  R  factor 
values  had  the  smallest  prediction  variances.  As  with  the  annual  R  factor  val¬ 
ues,  the  prediction  variances  were  higher  outside  of  the  rainfall  stations’  area. 
Within  the  rainfall  stations’  area,  the  sparser  the  rainfall  stations  and  the  higher 
the  spatial  variability  of  seasonal  R  factor  values,  the  larger  the  prediction  vari¬ 
ances.  These  features  for  the  large  area  could  also  be  applied  to  the  seasonal 
prediction  and  variance  maps  for  the  Fort  Hood  area  (Figure  8,  right).  The  esti¬ 
mates  and  their  variances  were  more  smoothing,  however,  because  Fort  Hood  is 
only  a  small  area  within  the  model  development  region. 

Figure  9  shows  the  24  prediction  maps  produced  from  the  spatial  and  temporal 
predictions  for  half-month  rainfall  erosivity  R  factor  values  for  the  large  area. 
The  predicted  R  factor  values  had  similar  spatial  and  temporal  distributions  to 
those  of  the  corresponding  observations.  Generally,  the  largest  values  were  at 
the  central  east  and  the  smallest  were  at  the  west.  The  spatial  variability,  how¬ 
ever,  was  more  smoothing  than  that  in  the  annual  and  seasonal  prediction  maps. 
The  predicted  half-month  R  factor  values  varied  from  0  to  85.  Overall,  the  val¬ 
ues  increased  slowly  over  the  first  six  half-months  and  rapidly  until  the  first  half 
of  June.  The  values  then  fluctuated  and  decreased  slightly  until  October;  then 
quickly  decreased  to  December. 

From  the  large  area  maps  in  Figure  9,  24  prediction  maps  of  half-month  rainfall 
R  factor  values  for  Fort  Hood  area  were  extracted  and  are  presented  in  Figure 
10.  The  predicted  values  varied  from  0  to  31.5  and  had  less  spatial  variability 
compared  with  those  for  the  whole  area.  They  rose  from  the  first  half  of  January 
to  the  first  half  of  June,  then  decreased  up  to  December.  Fluctuations  occurred 
in  some  months  (e.g.,  the  second  half  of  March  had  higher  values  than  the  first 
half  of  April,  and  the  second  half  of  December  had  larger  values  than  the  first 
half  of  December).  From  the  second  half  of  June  to  the  second  half  of  October, 
the  temporal  variability  of  the  predicted  values  was  relatively  small. 

The  prediction  variance  images  of  half-month  rainfall-runoff  erosivity  R  factor 
values  for  the  large  area  were  created  but  are  not  shown  here.  As  observed  in 
spatial  variability  of  the  original  observations,  the  variances  increased  overall 
from  January  to  June,  then  decreased  up  to  December  except  for  the  time  from 
the  second  half  of  September  and  the  first  half  of  October,  in  which  the  variances 
were  greater  than  those  in  August  and  the  first  half  of  September.  The  largest 
variances  occurred  in  the  first  half  of  June.  The  24  variance  maps  for  the  Fort 
Hood  area  (Figure  1 1)  were  extracted  from  the  corresponding  prediction  variance 
maps  for  the  large  area.  The  variances  changed  from  0  to  233.  When  the  pre¬ 
dicted  values  were  small  and  their  spatial  variabilities  low,  their  uncertainty  in 
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terms  of  variances  was  low.  When  the  predicted  values  were  large  and  their  spa¬ 
tial  variabilities  high,  their  uncertainty  in  terms  of  variances  was  high.  The  first 
half  of  January,  for  example,  had  the  smallest  variances  because  of  small  rainfall 
R  factor  values,  and  the  first  half  of  June  had  the  greatest  variances  because  of 
large  R  factor  values. 


Figure  9.  Prediction  images  of  rainfall  erosivity  R  factor  values,  on  a  half-month  basis  over  the 


large  area  where  the  rainfall  stations  are  located. 
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Figure  10.  Prediction  images  of  half-month  rainfall  erosivity  R  factor  values  in  the 
Fort  Hood  area. 
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Figure  11.  Prediction  variance  images  of  half-month  rainfall  erosivity  R  factor  values 
for  the  Fort  Hood  area. 


The  differences  between  the  predicted  and  observed  R  factor  values  for  annual 
rainfall  erosivity  were  calculated  and  plotted  against  observations  using  the  test 
data  of  29  rainfall  stations  (Figure  12).  The  observations  varied  from  190  to  465, 
and  the  differences  obtained  using  the  simulation  varied  from  -120  to  120.  Over¬ 
estimation  occurred  in  the  areas  with  lower  R  factor  values  and  underestimation 
in  the  areas  with  higher  R  factor  values.  According  to  coordinates  of  30  rainfall 
stations,  the  rainfall  erosivity  R  factor  values  were  obtained  by  the  traditional 
method  of  interpolating  the  isoerodent  map  and  also  comparing  the  observations 
shown  in  Figure  12.  The  differences  fell  in  the  range  of  0  to  -250.  The  underes¬ 
timation  occurred  for  all  the  rainfall  stations  and  was  much  worse  in  direction 
and  amount  than  the  occurrence  of  underestimation  and  overestimation  by  the 
simulation. 
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Figure  12.  Difference  between  estimate  and  observation  errors  versus  observed  annual  rainfall 
erosivity  R  factor  values  for  the  large  area  using  simulation  and  a  traditional  isoerodent  map. 


The  simulation  also  resulted  in  over-  and  underestimations  to  the  seasonal  pre¬ 
dictions  in  Figure  13  and  the  overestimation  was  slightly  more  significant  than 
underestimation.  The  range  of  relative  errors  (the  ratio  of  a  difference  to  an  ob¬ 
servation)  for  the  seasonal  R  factor  values  was  the  smallest  in  summer,  followed 
by  autumn,  spring,  and  winter,  although  the  absolute  errors  were  smaller  in 
winter  than  in  other  seasons. 

The  error  analysis  used  the  test  data  of  29  rainfall  stations  for  the  annual  and 
seasonal  rainfall  erosivity  R  factor  values  (Table  3).  Overall,  the  simulation  for 
annual  R  factor  values  led  to  the  average  estimate  falling  into  the  confidential 
interval  at  the  significant  level  of  5  percent,  although  it  was  slightly  overesti¬ 
mated.  The  traditional  method  of  using  an  isoerodent  map  resulted  in  serious 
underestimation  with  a  large  and  negative  mean  difference  and  an  average  es¬ 
timate  out  of  the  confidential  interval.  The  root  mean  square  error  (RMSE)  us¬ 
ing  the  isoerodent  map  was  251  percent  of  that  using  the  simulation.  In  Table  3, 
the  average  estimates  for  summer  and  spring  rainfall  R  factor  values  were 
within  the  confidence  intervals  at  the  significant  level  of  5  percent,  and  the  cor¬ 
responding  average  estimates  for  winter  and  autumn  were  very  close  to  the  up¬ 
per  limit  of  the  confidential  interval. 
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Figure  13.  Difference  between  estimates  and  observations  versus  observed  seasonal 
rainfall  erosivity  R  factor  values  for  the  large  area  using  the  simulation. 


Table  3.  Error  analysis  for  predicting  annual  and  seasonal  rainfall  erosivity  R  factor  values 
using  the  test  data  from  29  observations. 


R  factor 

Mean  of 
Test  Data 

Stdv  of 
Test  Data 

Confidential  Interval 

Mean  of 
Estimates 

Difference 

RMSE 

Lower 

Upper 

Annual 

359.06 

71.60 

331.60 

386.66 

373.34** 

14.28** 

212.17* 

-146.9* 

Winter 

30.04 

16.02 

23.87 

36.22 

37.37 

7.33 

Spring 

90.19 

28.64 

79.15 

101.23 

100.22 

10.02 

Summer 

127.22 

26.40 

117.04 

137.40 

125.12 

-2.10 

Autumn 

91.54 

28.42 

80.59 

102.50 

103.81 

12.27 

Note:  Stdv  =  standard  deviation,  RMSE  =  root  mean  square  error,  the  significant  level  is  5  percent, 
*  =  traditional  method,  and  **  =  simulation. 


The  estimates  and  observations  for  time  series  of  half-month  rainfall  R  factor 
values  were  compared  using  the  test  data  in  Figure  14.  All  the  average  esti¬ 
mates  fell  within  the  confidential  intervals  and  were  very  close  to  their  mean  ob¬ 
servations.  The  differences  between  the  estimated  and  observed  values  fluctu¬ 
ated  around  zero  from  -2.4  to  1.8.  The  maximum  RMSE  existed  in  the  first  half 
of  October  and  March,  when  the  largest  variation  of  the  rainfall  R  factor  was 
observed  in  the  test  data. 
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6  Conclusions 


The  sequential  Gaussian  simulation  strategy  reported  in  this  study  provided  the 
spatially  and  temporally  predicted  values  and  their  uncertainty  measures  in 
terms  of  prediction  variances  for  rainfall  erosivity  R  factor  values  in  predicting 
soil  loss  at  unknown  locations  and  areas.  The  spatial  and  temporal  distributions 
of  the  predicted  values  were  similar  to  the  observed  data  from  the  rainfall  sta¬ 
tions.  This  method  reproduced  the  variability  of  the  data  over  space  and  time, 
and  can  thus  be  recommended  as  a  monitoring  and  mapping  strategy  for  spatial 
and  temporal  prediction  and  uncertainty  analysis  of  the  rainfall  erosivity  R  fac¬ 
tor  in  predicting  soil  loss.  The  results  also  provide  important  inputs  required  to 
complete  an  overall  uncertainty  analysis  of  the  RUSLE. 

The  rainfall  erosivity  R  factor  is  an  important  variable  in  the  prediction  of  soil 
loss.  However,  it  is  difficult  to  derive  the  R  factor  values  in  areas  where  no  rain¬ 
fall  stations  exist.  Traditionally,  the  most  widely  used  method  is  to  interpolate 
the  R  factor  values  from  the  isoerodent  maps  where  R  factor  values  are  assumed 
constant  over  time.  The  method  presented  in  this  report  suggested  a  possible 
improvement  over  the  traditional  approach  in  deriving  the  R  factor  values. 

In  fact,  the  results  from  the  validation  data  showed  that  all  the  average  esti¬ 
mates  of  the  R  factor  by  the  simulation  for  annual,  seasons,  and  half-months 
were  within  the  confidential  intervals,  except  that  winter  and  autumn  data  were 
very  close  to  the  upper  limit.  The  average  estimate  of  annual  rainfall  R  factor 
values  by  the  isoerodent  map  was  outside  its  corresponding  interval  and  system¬ 
atically  negative.  The  annual  rainfall  R  factor  values  obtained  by  the  simulation 
for  the  Fort  Hood  areas  without  any  rainfall  stations  varied  from  350  to  376,  fal¬ 
ling  into  the  R  factor  values  of  four  rainfall  stations  around  it,  but  much  higher 
than  the  R  factor  of  270,  which  was  based  on  the  isoerodent  map.  That  is,  the 
annual  R  factor  values  from  the  traditional  isoerodent  map  for  the  large  and 
small  areas  were  much  less  than  the  observed  and  simulated  values,  which  is 
consistent  with  the  findings  by  McGregor,  Mutchler,  and  Bowie  (1980).  The  rea¬ 
son  for  the  lower  traditional  values  may  be,  in  part,  because  different  equations 
and  data  sets  were  used  for  calculation  of  the  R  factor  values  for  each  rainfall 
station. 
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The  results  in  this  study  also  implied  that  the  annual  rainfall  erosivity  R  factor 
had  a  large  variability  over  space.  Even  within  a  relatively  small  area,  such  as 
Fort  Hood  with  an  area  of  87,890  ha,  the  spatial  variability  may  not  be  over¬ 
looked.  Moreover,  there  was  a  high  temporal  variability  of  the  R  factor  in  the 
time  series  of  seasons  and  half-months.  As  expected,  the  summer  had  the  larg¬ 
est  R  factor  values,  then  autumn,  spring,  and  winter.  The  half-month  rainfall  R 
factor  values  increased  from  January  to  June,  then  fluctuated  and  decreased 
slowly  to  October,  and  after  that  tended  to  rapidly  decrease  to  December.  These 
findings  indicate  a  necessity  to  derive  the  local  estimates  and  their  uncertainties 
in  both  space  and  time  for  the  soil  loss  prediction  system  in  which  spatial  and 
temporal  variability  cannot  be  overlooked.  On  the  other  hand,  these  findings 
also  imply  an  importance  of  vegetation  cover  to  reduce  soil  loss  in  summer  by 
reducing  water  runoff. 

When  an  isoerodent  map  is  used  to  estimate  the  rainfall  R  factor  value,  its  un¬ 
certainty  is  unknown.  The  simulation  method  reported  here  provides  estimates 
with  variances  at  any  unknown  locations.  Where  the  rainfall  stations  used  for 
model  development  were  dense  and  the  rainfall  R  factor  value  was  low,  small 
variances  resulted;  otherwise  the  variances  were  large.  Decisionmakers  can 
therefore  apply  the  rainfall  R  factor  estimates  based  on  careful  assessment  of 
their  uncertainties. 
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